Relaxation of Surface Profiles by Evaporation Dynamics 
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We present simulations of the relaxation towards equilibrium of one dimensional steps and si- 
nusoidal grooves imprinted on a surface below its roughening transition. We use a generalization 
of the hypercube stacking model of Forrest and Tang, that allows for temperature dependent next- 
nearest-neighbor interactions. For the step geometry the results at T = agree well with the t 1 ^ 
prediction of continuum theory for the spreading of the step. In the case of periodic profiles we 
modify the mobility for the tips of the profile and find the approximate solution of the resulting free 
boundary problem to be in reasonable agreement with the T = simulations. 

PACS-numbers:68.10.Cr,05.70.Ln,68.35.Ct,68.35.Md 



I. INTRODUCTION 

The relaxation of corrugated crystal surfaces above 
their roughening transition is well understood in terms of 
continuum theory [jD. Experiments [pj and Monte Carlo 
simulations ||[| on the decay of sinusoidal grooves agree 
quantitatively with the predictions of continuum theory 
and even deviations due to the anisotropy of the surface 
free energy have been observed Q . For profiles imprinted 
on a crystal facet below its roughening transition things 
are less settled. Several predictions of continuum theories 
P-|io[ for time and wavelength dependence of the decay 
exist, where differences are mainly due to the different 
treatment of the singularity in the surface free energy 
that emerges for the orientation below the roughening 
temperature Tr. Experiments where surface diffu- 
sion prevails, show trapezoidal profiles with flat tops and 
bottoms, a feature qualitatively reproduced by several 
continuum approaches J8|,0. Simulations below Tr are 
severely hampered by the slow decay kinetics, a prob- 
lem somewhat less important for evaporation dynamics. 
Furthermore for T < Tr lattice effects play a role which 
are not present for T > Tr and which are not taken into 
account in the continuum theory. To suppress these ef- 
fects, larger systems need to be simulated. The scope of 
this article is to present new Monte Carlo simulations for 
the case of evaporation dynamics, which are able to test 
the predictions of continuum theory. In chapter II wc 
introduce a modified version of the hypercube stacking 
model |[l| for the case of a two dimensional surface. In 
chapter III we recall the results of continuum theory for a 
bunch of straight steps and compare them to simulations 
at T = 0. In chapter IV we show how one can modify the 
continuum theory for periodic grooves to get agreement 
with our T = simulations and discuss another recent 
attempt to this problem. In chapter V we give a con- 
clusion of our findings. Appendix A briefly describes the 
core part of our modifications of the hypercube stacking 



model. 



II. MODIFIED HYPERCUBE STACKING MODEL 



We use the hypercube stacking model of Forrest and 
Tang, which is a solid on solid model of the (111) surface 
of a simple cubic crystal as described in p] . Originally 
Blote and Hilhorst |l2| considered an anti-ferromagnetic 
Ising model on a triangular lattice with anisotropic cou- 
plings Ki = pE =: K, K 2 = K + (3E 2 and K 3 = 
K + {3Es in the limit K — > -co. They noted already that 
in two dimensions each Ising-spin configuration can be 
mapped on a SOS height configuration. In the isotropic 
case K\ = K 2 = K 3 considered in |ll]] each configu- 
ration is a ground state of an anti ferromagnetic Ising 
model with one frustrated bond per triangle and the cor- 
responding SOS model can be considered as temperature 
independent (E 2 = E$ = 0) or as being at infinite tem- 
perature (/3 — > but K — > — oo). All possible surface 
configurations have the same energy, since the number 
of broken nearest neighbor bonds is conserved. Dynam- 
ics can be imposed by assigning the same rate po to all 
evaporation and condensation events allowed by the SOS 
restriction. Spin flips which maintain the condition of 
having one frustrated bond per triangle correspond to 
adding or removing an atom without violating the SOS 
condition. 

We now introduce next-nearest-neighbor interactions of 
strength J. If an atom is removed from the surface 
the surface energy changes by AE = An = J(n — 3), 
where An is the difference in the number of broken next- 
nearest-neighbor bonds and n G {0, ...,6} is the number 
of next-nearest-neighbors in the same layer. Thus the to- 
tal surface energy is proportional to the number of bro- 
ken next-nearest-neighbor bonds. We impose the usual 
Metropolis rates 



f PoQ n ~ 3 if n > 3 
\ Pq if n < 3 



(1) 
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with q = exp(— -^y) for an evaporation event and 



w 



J p q 3 " if n < 3 
\ Po if n > 3 



(2) 



for a condensation event. Note that these rates fulfill de- 
tailed balance. 

As shown in detail in |TT| the spin representation of the 
model can be programmed very efficiently by using a 
multi-site coding algorithm, where each spin is coded by 
a single bit of an integer variable and the dynamics is 
incorporated via logical operations on the integers. In 
this type of algorithm we can also incorporate the next- 
nearest-neighbor interactions between surface atoms as 
will be shown in some detail in Appendix A. We use four 
dynamical sub-lattices as in |]l3| for the evolution of the 
height profile. The active sub-lattice is chosen at random 
and each evaporation or condensation event on that sub- 
lattice allowed by the SOS restriction takes place with 
the metropolis rate p of the event determined by eqs. (|l|) 
and (0) with po — \- Our generalization of the model 
of leads to a fast and memory saving implementa- 
tion of a solid on solid model, which we use in the sequel 
to study the decay of surface corrugations consisting of 
many steps. We should note that, due to the multi-site 
coding, where we use each bit of a random number to up- 
date one spin, the values of q, and therefore the available 
temperatures, are restricted to q 6 {0, k/2 rn , 1} with m 
integer and k odd integer between 1 and (2 m — 1) JlT| . 
The best performance of the algorithm is achieved for 
small m since 3m + 1 random integers are needed for 
each update of the spins stored in one integer. The cases 
q G {0, 1} (corresponding to T = and T = oc) need 
only one random number per update for the rate po = | 
while q = \ (corresponding to T = lAAJ/ks) needs al- 
ready four, namely one for Pq = \ and three for q 3 = g. 
After completion of this work we became aware of ref. 
p3f , where the same kind of model (interpreted as a fee 
crystal) with different rates is studied for periodic grooves 
at a temperature T — 0.68T# (corresponding to q = |). 
In the sequel we concentrate on the case T = 0. Then 
evaporation or condensation events can take place only 
at (11) steps since adding or removing an atom there 
does not increase the number of broken bonds. Decay 
kinetics is slowed down considerably by this choice, but 
as an advantage nucleation of islands or holes on ter- 
races is entirely absent. This simplifies the comparison to 
continuum theory, since now the mobility of the surface 
is proportional to the kink density. The (11) steps are 
rough even at T = due to the random update of growth 
and evaporation sites and the interaction between steps 
in this limit is purely entropic. We study the relaxation 
of two different initial profiles, a train of equidistant (11) 
steps separating two plane (111) surfaces (step geome- 
try) and a periodic sinusoidal corrugation of wavelength 
L. For the second we grow an initial profile that satisfies 



h(x,0) < aint(sin( (where aint denotes the integer 
part of the argument) on a flat surface of size L x L so 
that the steps of the profile do have (11) orientation. For 
the original model which is effectively at T = oo we re- 
cover the well known results of Mullins Q valid above 
the roughening transition. 



III. STEP GEOMETRY 

A. Continuum theory 

To describe the flattening of a profile under evapora- 
tion dynamics with a continuum theory we use the usual 
assumption that the profile evolution minimizes the free 
energy of the surface in the most direct way: 



a sf 
_ /l = _ M(u) _ 



(3) 



Here n(u) is the mobility of a surface of slope u = 1^ and 



dx 



the free energy F of the surface is given by the integral 
over the free energy a(u) per unit projected surface area. 
Below the roughening transition a(u) is written as uS 



a(u) — cr + <Ti\u\ H er 7 |u| 



(4) 



with 7 = 3 for purely entropic step repulsion. a§ de- 
scribes the energy of the (111) surface, a\ is the step free 
energy which contributes proportional to the step den- 
sity \u\ and the a 7 term describes step step interactions. 
Now the wedge singularity of the o\ term at u = leads 
to a (^-function singularity in (Q) . In the case of evapora- 
tion dynamics well below Tr, where thermally activated 
nucleation is negligible and the whole dynamics is due to 
step motion, the mobility of a surface is proportional to 
the step density j7p 



(5) 



Since nucleation is an activated process, it is entirely ab- 
sent in our T — simulations and eq. (JSJ) is valid even for 
a driven surface. Now the vanishing mobility at u = 
for the high symmetry orientation cancels the 5-function 
singularity in @ leading to 



5' 



/iio~(7 - l)|up 



dx 2 



■h . 



(6) 



For the step geometry M the scaling ansatz h(x,t) = 
Q(t~ a x) gives a — 1/(7 + 1) and leads to a Barenblatt 
solution |l5||l6|] for the slope 



(7) 



u(x,t) = c t- a (a 2 -y 2 ) 1 ^- 1) 



where y — x/t a , cq and a are constants determined by 
the initial condition and + denotes the positive part of 



2 



the bracket. Since eq. is invariant under the rescaling 
x = x/L, h = h/L and t — t/L 2 , the whole L dependence 
for a set of initial profiles of the same slope u but of dif- 
ferent size L can be absorbed by rescaling. This leads to 
cq ~ L°. and a 2 ~ £ 2 -4a j? or ^ ne coefficients of eq. (|7|). 



B. Simulational results 

How do these results compare with our T = simu- 
lations? Figure 1 displays the evolution of a step train 
of 11 steps on a 360 x 360 lattice up to t = 10 4 Monte 
Carlo steps with At = 500. The height is averaged over 
the columns parallel to the average step direction and 
over 10 3 independent runs. To check eq. (0) and the 
L dependence of the parameter cq and a we plot the 
scaled slope u ■ t 1 / 4 against the scaled width x ■ t -1 / 4 us- 
ing 7 = 3. Figure 2 displays the flattening data for three 
step trains of 45, 22 and 11 steps with the same initial 
slope. The topmost curves are data for 45 steps on a 
1440 x 1440 lattice up to t = 10 4 with At = 10 3 , aver- 
aged over 100 independent runs. The middle set gives the 
same information for 22 steps on a 720 x 720 lattice with 
At = 500 averaged over 10 3 runs and the lowest set are 
the differentiated data of figure 1. The curves of figure 
2 nicely approach the ellipses of the Barenblatt solution 
(fjj) and a good fit can be achieved with Co = 0.028 and 
a 2 = 1024, 500.4, 250.3 (from top to bottom) in agree- 
ment with Co ~ L and a 2 ~ L. The slight skewness in 
the initial condition (due to computational convenience) 
vanishes in course of time, indicating that the Barenblatt 
solution is attractive also for asymmetric initial profiles. 
The singularity of eq. (0) for u = is smeared out by the 
fluctuations of the rightmost and leftmost steps. Further- 
more one can resolve the mean positions of the individual 
steps, as is seen best for the 11 step data. Similar fea- 
tures have been observed recently jl7| for two steps with 
surface diffusion dynamics. 



IV. PERIODIC GROOVES 

A. Continuum theory 

For periodic surface grooves eq. (0) leads to a solution 
that predicts a non parabolic sharpening of the profile 
tips proportional to (5x)^ +1 ^ J J?],|| that was observed 
neither in previous p|,|l^-]n| nor in the present simula- 
tions. As pointed out by Rettori and Villain M the decay 
of the profile tips proceeds via the shrinking of islands 
formed by the two topmost meandering steps, which an- 
nihilate each other on contact. This process induces a 
nonzero mean kink density at the tips of the profile and 
hence the mobility does not vanish for u = 0. In this case 
we expect the wedge singularity of eq. (|^) to be rounded 



to an analytical function with nonzero /i(0) = /io- Dif- 
ferentiating eq. (|^) with respect to x and performing the 
functional derivate we get 



(8) 



dv .du. 2 , .d 2 v,du 2 dv d 2 u 

= /i(u) ^ ( fe } 

where v = a (u). Provided that fi(u) contains no L 
dependence eq. (g) is independent under the rescaling 
x = x/L, h = h/L and i — t/L 2 . Thus the data for a set 
of initial profiles with the same slope but different wave- 
lengths L should collapse on a single scaling curve under 
the rescaling. Note that the first two terms in are 
positive while the third one is negative. In the sequel we 
restrict the discussion to the interval x £ [0,L/4], where 
u > 0, which is sufficient to describe the whole profile by 
symmetry and periodicity. After initial transients have 
died out we expect dtu < for all u > 0, so the third term 
is dominant in eq. (|3) for u > 0. For u — > and 7 > 2 
the first term in eq. (B|) is much smaller than the second 
since /i'(tt) — > 0. In the step flow regime, where eq. (^) 
is valid and the third term dominates, the first term is 
of the same magnitude as the second. By neglecting the 
first term in eq. (0) we get 



d 

di' 1 



H{u) 



d 2 v 
dx 2 



(9) 



Now the singularity in the free energy is not canceled by 
a vanishing mobility and gives rise to a 5' function in ^3 . 
Bonzel and Preuss jl0| smoothened the wedge singular- 
ity and solved the resulting eq. numerically in the case of 
surface diffusion. This leads to qualitatively correct re- 
sults but introduces additional parameters. As required 
by the thermodynamic stability of neighboring surface 
orientations with a / 0, d = fJ '(' u ) is strictly increasing 
and therefore invertible. So we can transform the singu- 
lar eq. Q) into a well defined free boundary problem § 
by using v instead of u as the independent variable. Note 
that for a singular orientation with unstable neighbor ori- 
entations |2Q| such a transformation cannot be applied. 
Inverting v = u'{u) we obtain || 

^1/(7-1) if v <-a 1} 
u(v) = I if - 01 < v < ai , (10) 
t (ifc-trO) 1 /^-!) if v>a x . 

By inserting ( |To| ) into (^) we arrive at the free boundary 
problem: 

d 1 7 _ 2 Q 2 CL 

—v = (j - l)^)^- 1 (v - a!)^ -^v for ig [0, — ] 

(11) 
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and 



B. Simulational results and discussion 



dx 2 



= for 



(12) 



where C £ [0, 1] is a time dependent free boundary. Here 
C < 1 means the existence of facets of size 2(1 — C,)L at 



the tips of the profile. The solution of eq. (12) taking 
into account point symmetry is 



v(x, t) = G\ 



L-Ax 



L-((t)L ■ 
The boundary conditions for eq. ( |TT| ) are 
4cti 



(13) 



v{ -> t) = -L-C(t)L 



and v'(0,t)=0. (14) 



Eq. dTJ) cannot be simply solved by a scaling ansatz as 
eq. (|6j) since this ansatz does not fulfill the boundary con- 
ditions. Furthermore we don't know the precise form of 
the mobility fi(u) near the tips. Thus we have to try an 
approximative or a numerical solution. In the step flow 
dominated region of not too small slope, which for our 
T = simulations covers the whole profile except the tips 
where step step annihilation takes place we can safely use 
cq. (||) for the mobility. At the facet edges however the 
solution is determined by the boundary conditions. As 
a first approximation we use a power series ansatz up to 
second order in x, namely v — <J\ = ao(t)(l — ^ffr), that 
fulfills the boundary condition with ao(t) — <7i£/(2 — 2£) 
but solves eq. ([ll]) up to second order in x only in the 
singular limit £ — > 1. For ( « 1 we approximately find 
'Mi (7 — l)ci/£ 2 = : — ci, which leads to 



C = c - cit , 



(15) 



with a constant cq determined by the initial condition. 
For the amplitude of the profile we find 



h(—,t) = — (— i-)^ 1 ( C °- Cl ^ 



4 "2a 



7 (1 — cq + cit) t- 1 Jo 



(l-y 2 y dy 
(16) 



and the profile shape for 7 = 3 is 



h{x, t) 



2, XL ..Ax 
7r 4 QL 



(t^) 2 + arcsin{-^)) 



(17) 



We note that the solution up to this order is of scaling 
form and the time dependence enters only via By 
adding a term ai(t){\ — ^Ify) 2 we can solve eq. ( ]TT| ) up to 
fourth order in x for arbitrary If we assume <C a 
we can solve the emerging differential equation for £ nu- 
merically which gives the time dependence of the decay 
and fulfills the assumption a\ ao self-consistently. For 
(«lwe recover eqs. (p3|)-(|i"7|). 



Figure 3 displays simulational data of initial sinusoidal 
grooves on square lattices. The corresponding values of 
system size L, initial height ho, number of Monte-Carlo 
steps t and number of independent runs are given in table 
1. The data have been averaged in the transverse direc- 
tion and over the independent runs. We fitted the pro- 
files with eq. ([l7]) taking the actual amplitude as Mff) 
and used ( as fit parameter. The values for the abso- 
lute deviations Ahfu of the fits and for Qn are also dis- 
played in table 1. We can not decide wether the remain- 
ing small, but systematical deviations from eq. ( J17|) are 
due to higher order terms of the solution of eq. fllip or 
caused by our approximation leading to eq. (^). As for 
the step trains the singularity at the facet edge and also 
the facet itself are blurred by step fluctuations. Since £ 
tends to 1 for larger L wc don't expect to find a macro- 
scopic facet for L — > 00. Eq. (|l^) which describes the 
time dependence of the amplitude could not be tested 
seriously with the data of figure 3 since the time evolu- 
tion was too short to pin down the three independent 
parameters of eq. ( ^6|) p recisely. Nevertheless the values 
of Co and c\ via eq. (|15[) give an alternative estimate for 
the time evolution of £ which shows a considerably faster 
decay than the direct evaluation via eq. (|l7|). This in- 
dicates that eqs. (p^|)-(pr|) don't tell the full story and 
higher order terms in the solution of ([ll]) are necessary 
to describe the full time dependence of the decay. 
To test the wavelength dependence of the decay we in 
figure 4 plot the scaled amplitude h — h/L against the 
scaled time t — t/L 2 . The upper 5 curves of figure 4 
display the scaled amplitude for L = 1440, 720, 540, 360, 
180 from top to bottom. The lower set of curves shows 
the same sequence of data with a different time scaling 
t/L 2 - 15 . The obvious violation of the expected scaling 
cannot be traced back to our approximations, since al- 
ready eq. (||) is invariant under reparametrisation. Wc 
can think of at least two possible sources. We start with 
a set of straight steps which roughen at the begin of the 
decay. This feature is not present in our continuum ap- 
proach and might introduce an additional L dependence. 
Secondly the step step annihilation at the tips of the 
profile may induce a L dependence of the mobility not 
present in the case of a step train. Actually Tang in his 
recent work did find evidence for logarithmic corrections 
to the expected scaling due to the top step annihilation 



process |23 



To describe the profile form Tang modifies the equa- 
tion of Langon and Villain by adding a term pro- 
portional where the function so(t) is determined 
self-consistently out of the width of the topmost terrace. 
This approach, leading to fits for the profile shape of 
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an accuracy comparable to ours, effectively introduces a 
mobility /i(0) / as we did, but in contrast to our de- 
scription Tang ignores the wedge singularity of the free 
energy which led us to the free boundary problem. Thus 
Tang ^jl conjectures that the equation of Lanqon and 
Villain becomes exact in the limit L — > oo. In contrast 
we expect the solution of the free boundary problem i.e. 
cqs. (|l5|)-(|l7|) plus small corrections due to higher orders 
to be valid for L — > oo. The different predictions for the 
profile shape in the limit L — > oo are displayed in figure 
5. One finds that the theory of Langon and Villain pre- 
dicts more pronounced tips than the solution of the free 
boundary problem for both 7 = 2 and 7 = 3. Our largest 
simulations of wavelength 1440 support eq. (O) but they 
presumably can also be fitted with the term induced by 
Tang Q. Further effort is necessary to establish the 
wavelength dependence for periodic grooves. 



V. CONCLUSION 

In conclusion we showed that our large scale T = 
simulations of step trains agree quantitatively with the 
predictions of continuum theory below the roughening 
transition. For periodic grooves we used a nonzero mo- 
bility /u(0) in the continuum theory to capture some of 
the subtleties of step step annihilation. Our approximate 
solution for the profile shape is in good agreement with 
the simulations while for time decay and wavelength scal- 
ing discrepancies remain. 

One can extend our simulation to temperatures T > to 
study the influence of nucleation on surface free energy 
and mobility. One can also incorporate surface diffusion 
dynamics in our kind of algorithm. This would be in- 
teresting since the existing simulations I H^]^ use too 
small systems for a comparison to continuum theory. 
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APPENDIX A: 

Here we give the logical function that counts the num- 
ber of atoms on the six in plane next-nearest-neighbor 
sites for a given evaporation or condensation site. This 
function is the main extension of the code given already 
in Jn]]. Let is(i) G {0,1} denote the value of the spin 
on an evaporation or condensation site and is(j) with 



j £ {1, ...,6} in some arbitrary order, the values of the 
spins on the six next-nearest-neighbor sites of site i. We 
define the flag fj :— is(j).xai.is(i), which has the values 



fj 



if h(j) — h(i) (occupied site j) 

1 if h(j) = h(i) — 3 (empty site j) 



for i being an (occupied) evaporation site, and 
if h(j) — h(i) (empty site j) 



1 if h(j) = h(i) + 3 (occupied site j) 



(Al) 



(A2) 



for i being a (empty) condensation site. Thus fj is the 
occupation number of site j if site i is a condensation 
site and .not. fj if i is a evaporation site. Out of the fj 
we now create flags nu with the value 1 for k occupied 
next-nearest-neighbor sites, being otherwise: 

n = .not.(/i.or./ 2 .or./3.or./4.or./ 5 .or./ 6 ), (A3) 

rii = (/i.xor./ 2 .xor./ 3 .xor./ 4 .xor./ 5 .xor./ 6 ).and.not. 
{[(/i.and./ 2 ).or.(/3.and./ 4 ).or.(/5.and./ 6 )] 



i[Ui-and./ 2 ).or.(/3.and./ 4 ).or.(/5.and./ 6 )J .or 
[(/i.and./ 2 ).or.(/3.and./5).or.(/4.and./ 6 )] .or. 
[(/i.and./ 2 ).or.(/ 3 .and./ 6 ).or.(/4.and./ 5 )]}, 



(A4) 



ni = (/i.or./ 2 .or./ 3 .or./ 4 .or./ 5 .or./ 6 ).and. [.not. 

(/i.xor./ 2 .xor./ 3 .xor./ 4 .xor./ 5 .xor./ 6 )] .and. [.not. 
{[(/i.and./ 2 ).or.(/ 3 .and./ 4 ).or.(/ 5 .and./ 6 )] .and. 
[(/ 1 .and./ 3 ).or.(/ 2 .and./ 5 ).or.(/ 4 .and./ 6 )]}] , (A5) 

71.4 = [.not.(/i.and./ 2 .and./ 3 .and./ 4 .and./5.and./e)] .and. 
[.not.(/i.xor./ 2 .xor./ 3 .xor./ 4 .xor./ 5 .xor./ 6 )] .and. 
{[(/i-or./ 2 ).and.(/ 3 .or./ 4 ).and.(/ 5 .or./ 6 )] .or. 
[(/ 1 .or./ 3 ).and.(/ 2 .or./ 5 ).and.(/ 4 .or./ 6 )]}, (A6) 

n 5 = (/i.xor./ 2 .xor./ 3 .xor./ 4 .xor./ 5 .xor./ 6 ).and. 

{[(/i-Or./ 2 ).and.(/ 3 .or./ 4 ).and.(/ 5 .or./ 6 )] .and. 
[(/i-or./ 2 ).and.(/ 3 .or./ 5 ).and.(/ 4 .or./ 6 )] .and. 
[(/ 1 .or./ 2 ).and.(/ 3 .or./ 6 ).and.(/ 4 .or./ 5 )]}, (A7) 



n e = (/i.and./ 2 .and./ 3 .and./ 4 .and./ 5 .and./ 6 ). (A8) 

Note that we listed the flags for the case of i be- 
ing a condensation site, and one has to replace all fj by 
.not./j in the case of i being a evaporation site. With 
these flags at hand we can construct expressions for the 
missing rt 3 and for n{ 3 6} = -not.(no.or.ni.or.n 2 ) and 
n {o....,3} = -not. {uq. or. n§. or. 714), which we need to apply 
the Metropolis rates (Q) and (||). 
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FIG. 1. Flattening of a step train of 11 initially equidistant 
steps on a 360 x 360 lattice up to t — 10 4 Monte Carlo steps 
with At = 500. 

FIG. 2. Time evolution of the scaled slope u ■ £ 1//4 plotted 
versus the scaled width x ■ i~ 1//4 for step trains of 45, 22, 11 
steps (from top to bottom) with the same initial slope. 

FIG. 3. Flattening sequences of periodic sinusoidal grooves 
with the same initial slope on quadratic lattices of L = 180, 
360, 540, 720, 1440. 

FIG. 4. The upper 5 curves display the time evolution 
of the scaled amplitudes h/L for the data of figure 3 plot- 
ted versus the scaled time t/L 2 . The lower set displays the 
same data plotted versus a differently scaled time t/L 2 ' 15 . 
In each set from top to bottom data for the system sizes 
L = 1440, 720, 540, 360, 180 are displayed. 

FIG. 5. Different predictions for the profile shape in the 
limit L — > co are plotted with the same initial slope at 
x = O.From top to bottom the curves are: numerical solu- 
tion of eq. (4) (after a separation ansatz h(x,t) = g(t)f(x)) 
with 7 = 3 (topmost) and 7 = 2, approximate solutions of 
the free boundary problem (with f = 1) for 7 = 3 (eq. (15)) 
and 7 = 2, and a sine-function for comparison. 



TABLE I. Parameter for the simulational data of figure 3 



L 


h 


t x 10 4 


runs 


Cfit 


Aft/it 


180 


11.85 


1.5 


10 3 


0.89-0.93 


0.05 


360 


23.7 


4 


200 


0.96 


0.1 


540 


35.5 


4 


100 


0.97 


0.1 


720 


47.5 


10 


20 


0.98 


0.2 


1440 


94.5 


40 


5 


0.99 


0.4 



G 



0.070 
0.065 
0.060 
0.055 
0.050 
0.045 
0.040 
0.035 
0.030 
0.025 
0.020 




0.0000 



1 

0.0001 



+ 



+ 



0.0002 0.0003 
scaled time 



1 

0.0004 



1 

0.0005 



1.2 



1 



. 8 



. 6 



. 4 



0.2 



